

function f = fit_fundamental( matches,method )
%method true: eight-point algorithm
%method false: normalized eight-point algorithm
num = size(matches, 1);
pts1h=matches(:,1:2);
pts2h=matches(:,3:4);
if method
    m = zeros(num, 9);
    for idx = 1: num
      m(idx,:) = [...
        pts1h(idx,1)*pts2h(idx,1), pts1h(idx,2)*pts2h(idx,1), pts2h(idx,1),...
        pts1h(idx,1)*pts2h(idx,2), pts1h(idx,2)*pts2h(idx,2), pts2h(idx,2),...
                     pts1h(idx,1),              pts1h(idx,2), 1];
    end
    [~, ~, vm] = svd(m, 0);
    f = reshape(vm(:, end), 3, 3)';

    % Enforce rank-2 constraint
    [u, s, v] = svd(f);
    s(end) = 0;
    f = u * s * v';
    f = f / norm(f);
    if f(end) < 0
      f = -f;
    end
else
    % Normalize the points
    t1=get_normalization_matrix(pts1h);
    t2=get_normalization_matrix(pts2h);
    pts1h(:,3)=1;
    pts1h=pts1h*t1;
    pts2h(:,3)=1;
    pts2h=pts2h*t2;

    % Compute the constraint matrix
    m = zeros(num, 9);
    for idx = 1: num
      m(idx,:) = [...
        pts1h(idx,1)*pts2h(idx,1), pts1h(idx,2)*pts2h(idx,1), pts2h(idx,1),...
        pts1h(idx,1)*pts2h(idx,2), pts1h(idx,2)*pts2h(idx,2), pts2h(idx,2),...
                     pts1h(idx,1),              pts1h(idx,2), 1];
    end

    % Find out the eigen-vector corresponding to the smallest eigen-value.
    [~, ~, vm] = svd(m, 0);
    f = reshape(vm(:, end), 3, 3)';

    % Enforce rank-2 constraint
    [u, s, v] = svd(f);
    s(end) = 0;
    f = u * s * v';

    % Transform the fundamental matrix back to its original scale.
    f = t2' * f * t1;

    % Normalize the fundamental matrix.
    f = f / norm(f);
    if f(end) < 0
      f = -f;
    end
end


end


